df <- read.csv2(file = "data/ZNOandVoating/input_2017.csv", header = FALSE, sep = ",", dec = ".")
names(df) <- c("ukr","math", "pro_ukr", "radical", "opposition", "small", "not_voted")

head(df, 3)

##   ukr math  pro_ukr  radical opposition    small not_voted
## 1 110  128 48.82355 8.705518   0.344142 5.856787     36.27
## 2 119  124 48.82355 8.705518   0.344142 5.856787     36.27
## 3   0    0 48.82355 8.705518   0.344142 5.856787     36.27

ЗНО з математики та української мови

Виведемо загальні характеристики даних ЗНО з математики та української мови:

describe(df[, c("math", "ukr")])

## df[, c("math", "ukr")] 
## 
##  2  Variables      105801  Observations
## ---------------------------------------------------------------------------
## math 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       82        1    140.5    42.21        0      105 
##      .25      .50      .75      .90      .95 
##      121      145      169      185      189 
## 
## lowest :   0.0 100.0 101.0 102.0 104.0, highest: 197.5 198.0 199.0 199.5 200.0
## ---------------------------------------------------------------------------
## ukr 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       53    0.995    117.6    60.25        0        0 
##      .25      .50      .75      .90      .95 
##      107      129      157      177      185 
## 
## lowest :   0 100 103 107 110, highest: 196 197 198 199 200
## ---------------------------------------------------------------------------

Із таблиці вище видно, що в записах нема пропущених значень. Судячи із значень квантилів, незначна частка результатів ЗНО з математики є нулями (менше 10%). Варто відмітити, що трішки більша частка значень 0 у ЗНО з української мови (бфльше 10%, але менше 25%).

Подивимося на розподіли:

p1 <- plot_ly(x = df[, "math"], type = "histogram", name = "math")
p2 <- plot_ly(x = df[, "ukr"], type = "histogram", name = "ukr") 

plotly::subplot(p1, p2) %>% layout(title = "EIT: Ukranian language and mathematics")

remove(p1); remove(p2)

Розподіли балів результатів ЗНО з української та математики помітно різняться.

Подивимося на співвідношення людей, котрі отримали 0 з математики та українській одночасно.

print("Не подолали поріг або не здавали ЗНО відносно всіх абітурієнтів: " %&% as.character(round(nrow(df[(df$math == 0 & df$ukr == 0), ]) * 100 / nrow(df), 2)) %&% "%")

## [1] "Не подолали поріг або не здавали ЗНО відносно всіх абітурієнтів: 3.56%"

print("Не подолали поріг або не здавали ЗНО відносно тих, хто не склав математику: " %&% as.character(round(nrow(df[(df$math == 0 & df$ukr == 0), ]) * 100 / nrow(df[(df$math == 0), ]), 2)) %&% "%")

## [1] "Не подолали поріг або не здавали ЗНО відносно тих, хто не склав математику: 70.93%"

print("Не подолали поріг або не здавали ЗНО відносно тих, хто не склав українську: " %&% as.character(round(nrow(df[(df$math == 0 & df$ukr == 0), ]) * 100 / nrow(df[(df$ukr == 0), ]), 2)) %&% "%")

## [1] "Не подолали поріг або не здавали ЗНО відносно тих, хто не склав українську: 21.74%"

Із виведеної аналітики напрошуються наступні висновки, 1) Невелика частка абітурієнтів, що мають 0 одночасно з двох предметів. 2) Абітурієнти, котрі не склали іспит із математики, скоріш за все не склали й українську.

Поглянемо на кореляції між результатами ЗНО, політичними силами:

c <- cor(df, method = "pearson")
plot_ly(x = rownames(c), y = colnames(c), z = c, colors = "Greys", type = "heatmap") %>% 
  layout(title = "Pearson Correlation")

remove(c)

Відмітимо майже 0 кореляцію між результатами ЗНО з української мови та про-українським політичним напрямком та слабку кореляцію між результатми математики та української.

Видалимо стрічки, що відповідають спостереженням із 0 хочаб в одному ЗНО та подивимося на розподіл.

df <- df[(df$math != 0) & (df$ukr != 0) ,]
p1 <- plot_ly(x = df[, "math"], type = "histogram", name = "math")
p2 <- plot_ly(x = df[, "ukr"], type = "histogram", name = "ukr")

plotly::subplot(p1, p2) %>% layout(title = "Distribution of the passed EIT",
  annotations = list(list(x = 0.2, y = 1.05, text = "math", showarrow = F, xref = 'paper', yref = 'paper'),
                     list(x = 0.8, y = 1.05, text = "ukr", showarrow = F, xref = 'paper', yref = 'paper')))

remove(p1); remove(p2);

Відмітимо, що наразі частка людей, котрі погано здали українську мову зменшилася близько вдвічі (було в околі 900, стало в околі 400).

Гістограма розсіювання, в свою чергу, виглядає наступним чином:

plot_ly(data = df, x = ~math, y = ~ukr, mode = "markers", type = "scatter",  marker=list(opacity = 0.7)) %>% layout(title = "Passed EIT results")

Варто відмітити, що для абітурієнтів, що склали математику краще ніж 190 балів, помітна тенденція, як правило, кращого складання ЗНО з української.

Подивимося на розподіли по кожній політичній силі окремо (якщо конкретно обрана політична сила має для даного спостереження найбільший вплив, то дане спостереження відносимо до ціє політичної сили)

• математика

W <- as.data.frame(df[, -c(1, 2)])
W[, "max_value"] <- apply(W, 1, max)

all_plots <- list()
annotations_list <- list()
names_ <- names(W[, -ncol(W)])

for (i in seq_along(names_)) {
  name <- names_[i]
  d <- df[W[, name] == W["max_value"], ]
  if (nrow(d) != 0) {
    all_plots[[name]] <- plot_ly(x = d[, "math"], type = "histogram", name = name)
    annotations_list[[i]] <- list(x = 0.15*i, y = 1.005, text = name, showarrow = F, xref = 'paper', yref = 'paper')
  }
}

plotly::subplot(all_plots) %>% 
  layout(title = "Mathematics EIT according to political strength", annotations = annotations_list)

remove(d); remove(all_plots); remove(annotations_list); remove(i); remove(name)

• українська мова

all_plots <- list()
annotations_list <- list()

for (i in seq_along(names_)) {
  name <- names_[i]
  d <- df[W[, name] == W["max_value"], ]
  if (nrow(d) != 0) {
    all_plots[[name]] <- plot_ly(x = d[, "ukr"], type = "histogram", name = name) 
    annotations_list[[i]] <- list(x = 0.15*i, y = 1.005, text = name, showarrow = F, xref = 'paper', yref = 'paper')
  }
  
}

plotly::subplot(all_plots[["pro_ukr"]], all_plots[["not_voted"]]) %>% 
  layout(title = "Ukr EIT according to political strength", annotations = annotations_list)

remove(d); remove(all_plots); remove(annotations_list); remove(names_); remove(name); remove(i)

Відмітимо, що розподіл результатів ЗНО з математики трішки різниться длярегіонів з різними політичними вподобаннями. У той час, як розподіли результатів з української, якщо порівнювати відносно йомірних груп, до яких віднесли відповідні області абітурієнтів, мають дуже близьку природу.

for (name in names(W[, -ncol(W)])) {
  d <- df[W[, name] == W["max_value"], ]
  if (nrow(d) != 0) {WVPlots::ScatterHist(d, "math", "ukr", name, smoothmethod = "none", estimate_sig = FALSE, minimal_labels = TRUE)}
}

remove(d); remove(W); remove(name); remove(df)

Із останньої пари зображень здається, для регіонів проукраїнського напрямку менший розкид даних і таким чином більш лінійна залежність. Тим часо для регіонів, які відповідають регіонам із неголосуванням, більш помітна тенденція: більший бал з математики - кращий результат з української.

Результати голосування

Виведемо загальні характеристики кожної політичної сили:

df <- read.csv2(file = "data/ZNOandVoating/input_2017.csv", header = FALSE, sep = ",", dec = ".")
names(df) <- c("ukr","math", "pro_ukr", "radical", "opposition", "small", "not_voted")
describe(df[, -c(1, 2)])

## df[, -c(1, 2)] 
## 
##  5  Variables      105801  Observations
## ---------------------------------------------------------------------------
## pro_ukr 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       25    0.997    30.98    14.68    9.772   15.712 
##      .25      .50      .75      .90      .95 
##   19.371   33.749   41.339   48.824   53.872 
## 
## lowest :  9.114851  9.771840 15.712444 16.226912 19.371495
## highest: 43.008240 44.778925 48.823553 50.417952 53.872000
## ---------------------------------------------------------------------------
## radical 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       25    0.997    6.228    2.842    2.245    3.047 
##      .25      .50      .75      .90      .95 
##    3.834    5.944    8.414    8.912   10.428 
## 
## lowest :  1.726920  2.245021  3.046992  3.831828  3.834072
## highest:  8.886996  8.911950 10.084956 10.427880 11.555155
## ---------------------------------------------------------------------------
## opposition 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       25    0.997    5.261    5.231    0.497    0.497 
##      .25      .50      .75      .90      .95 
##    1.343    2.535   10.467   12.536   14.575 
## 
## lowest :  0.344142  0.416508  0.497000  0.953295  1.042825
## highest: 10.466742 11.615622 12.027133 12.535560 14.574912
## ---------------------------------------------------------------------------
## small 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       25    0.997    9.467    2.086    7.189    7.361 
##      .25      .50      .75      .90      .95 
##    8.366    9.167   10.890   13.113   13.113 
## 
## lowest :  5.856787  6.778902  7.189344  7.360584  7.496422
## highest: 10.890088 11.198572 11.346280 13.112736 13.519935
## ---------------------------------------------------------------------------
## not_voted 
##        n  missing distinct     Info     Mean      Gmd      .05      .10 
##   105801        0       25    0.997    48.06    10.99    30.00    35.15 
##      .25      .50      .75      .90      .95 
##    42.76    45.47    54.68    60.48    67.13 
## 
## lowest : 30.00 31.72 35.15 36.27 39.79, highest: 57.20 58.64 60.48 67.13 67.60
## ---------------------------------------------------------------------------

Одразу видно, що 25 унікальних значення в кожному стовпчику відповідають 24 областям України і окремо м. Київ. Варто відмітити, що відсутні області, із переважними радикальними або опозиційними силами. Також відмітимо що немає регіонів, у яких голоса віддані виключно за малі політичні партії.

Розподіл абітурієнтів по областях:

ab <- plyr::ddply(.data = df, .variables = .(pro_ukr, radical, opposition, small, not_voted), .fun = nrow)
plot_ly(ab, y = ~V1, type = "bar") %>% 
  layout(title = "Abiturients distribution according regions", 
         xaxis = list(title = "Region ID"), yaxis = list(title = ""))

Грубо кажучи, половина областей має приблизно в 2 рази більше абітурієнтів, ніж інша частина регіонів.

Із розподілу відсотків підтримки кожної політичної сили відповідно до регіону видно, що більшість регіонів є політично неактивними. Також варто відмітити, що розподіли голосів малих партії, у більшості випадків, збігається з розподілом голосів радикальної партії.

plot_ly(y = ab[, "pro_ukr"], type = "bar", name = "pro_ukr") %>%
  add_trace(y = ab[, "not_voted"], name = "not_voted") %>%
  add_trace(y = ab[, "radical"], name = "radical") %>%
  add_trace(y = ab[, "opposition"], name = "opposition") %>%
  add_trace(y = ab[, "small"], name = "small")

remove(ab); remove(df)

Визначення параметрів згладжування

Скористаємося 5-ти фолдовою кросс-валідацією для знаходження оптимальних параметрів згладжування для кожної моделі аналогічно до того, як це було зроблено для результатів ЗНО 2016 року.

# h_range <- c(seq(0.1, 1, 0.2), seq(1, 5, 0.5), 8, 10)
# 
# for (i in 1:5) {
#   print("***** " %&% as.character(i) %&% " ******")
#   set.seed(i)
#   cv_split <- cross_validation_split(train, k = 5)
#   GNW_cv_results <- GNWcv_across_h(h_range = h_range,
#                                  cv_df_split = cv_split,
#                                  X_colname = "math", Y_colname = "ukr",
#                                  W_colname = c("pro_ukr", "radical", "opposition", "small", "not_voted"),
#                                  use_parallel = TRUE)
#   write.csv(GNW_cv_results, "data/computation_results/2017cc_data/5cv_mean_std_" %&% as.character(i) %&% ".csv")
# }
# 
# # h_opt <- optimal_h(GNW_cv_results)
# remove(h_range); remove(i); remove(cv_split); remove(GNW_cv_results)

Результат середніх та диспресії збережемно в змінній \(GNW\_cv\_results\), а оптимальні \(h\), серед запропонованих для перебору, у змінну \(h\_opt\).

data_path <- "data/computation_results/2017cc_data/"
file_names <- dir(data_path) #where you have your files
file_names <- data_path %&% file_names

cv_results <- lapply(file_names, 
                     function(x){
                       df <- as.matrix(read.csv(x, row.names = 1))
                       l <- list(df=df , h=optimal_h(df))
                       return(l)
                       })

cv_results_h <- lapply(cv_results, function(x){x$h})
cv_results_mean_std <- lapply(cv_results, function(x){x$df})

print(cv_results_h)

## [[1]]
##         [,1]     [,2]     [,3]     [,4]     [,5]
## [1,] 506.736 8698.633 546.4858 22349.13 285.8498
## [2,]   8.000    0.700  10.0000     3.50   2.0000
## 
## [[2]]
##          [,1]     [,2]     [,3]    [,4]     [,5]
## [1,] 2391.061 36638.74 484.1853 14307.5 515.5184
## [2,]    2.000     8.00   8.0000     2.0   2.0000
## 
## [[3]]
##          [,1]     [,2]     [,3]     [,4]     [,5]
## [1,] 2315.868 13828.72 395.3039 16448.56 1462.592
## [2,]   10.000     2.00  10.0000     1.50    2.500
## 
## [[4]]
##         [,1]     [,2]     [,3]     [,4]     [,5]
## [1,] 917.807 5505.784 592.3602 39636.91 487.1869
## [2,]   1.500    8.000   8.0000     1.50   4.5000
## 
## [[5]]
##          [,1]     [,2]     [,3]     [,4]     [,5]
## [1,] 810.9435 17021.33 504.6474 10613.12 239.7613
## [2,]   5.0000     2.00   3.5000     5.00   1.5000

Порівняння результатів кросс-валідації та прогнозування на відкладеній вибірці

Порівняємо тепер результати кросс-валідації з результатми прогнозування на тестовій частині даних.

# train <- read.csv("data/computation_results/train2017.csv", row.names = 1, stringsAsFactors = FALSE)
# test <- read.csv("data/computation_results/test2017.csv", row.names = 1, stringsAsFactors = FALSE)
# h_opt <- cv_results_h[[1]][2, ]
# 
# GNW <- GeneralisedNadarayaWatson$new()
# GNW$train(X_train = as.numeric(train[, "math"]),
#           Y_train = as.numeric(train[, "ukr"]),
#           W_train = as.matrix(train[, -c(1, 2)]))
# GNW$run_cluster()
# prediction_with_acoeff <- GNW$predict_in_parallel(X_test = as.numeric(test[, "math"]),
#                                                   W_test = as.matrix(test[, -c(1, 2)]), h = h_opt)
# GNW$stop_cluster()
# 
# prediction <- prediction_with_acoeff[["prediction"]]
# A_coeff <- prediction_with_acoeff[["A_test"]]
# 
# write.csv(prediction, "data/computation_results/prediction2017v.csv")
# write.csv(A_coeff, "data/computation_results/A_coeff2017v.csv")
# 
# remove(GNW); remove(train)

Порахуємо зважену суму МНК для прогнозу зробленого на тестових даних.

train <- read.csv("data/computation_results/train2017.csv", row.names = 1, stringsAsFactors = FALSE)
test <- read.csv("data/computation_results/test2017.csv", row.names = 1, stringsAsFactors = FALSE)

prediction <- as.matrix(read.csv("data/computation_results/prediction2017v.csv", row.names = 1, stringsAsFactors = FALSE))
A_coeff <- as.matrix(read.csv("data/computation_results/A_coeff2017v.csv", row.names = 1, stringsAsFactors = FALSE))

WMSE <- weighted_MSE(Y_true = as.numeric(test[, "ukr"]), Y_predicted = prediction, A_coeff = A_coeff)
names(WMSE) <- rep("comp.") %&% as.character(1:5)
WMSE

##        comp.1        comp.2        comp.3        comp.4        comp.5 
## -2.329684e+03 -2.176772e+05  1.076314e+02 -5.490781e+10  3.212800e+02

Відмітимо, що найменше значення зваженого МНК для моделей, що відповідають опозиційним силам.

Порівняємо зважені МНК результатів прогнозування з результатами крос-валідації, використовуючи відповідними параметрами згладжування.

h_opt <- cv_results_h[[1]][2, ]
GNW_cv_results <-cv_results_mean_std[[1]]

cv_mean_results <- sapply(1:length(h_opt),
                          function(i){
                            GNW_cv_results[which("mean_" %&% as.character(h_opt[i]) == rownames(GNW_cv_results)), i]
                            })
cv_test_comparing <- rbind(cv_mean_results, WMSE)
rownames(cv_test_comparing) <- c("cross-validation", "test")
colnames(cv_test_comparing) <- 1:5
print(cv_test_comparing)

##                          1           2        3             4        5
## cross-validation   506.736    8698.633 546.4858  2.234913e+04 285.8498
## test             -2329.684 -217677.232 107.6314 -5.490781e+10 321.2800

Із останньої таблиці видно, що квадрати залишків 3-ї та 4-ї моделей, порахованих для кросс-валідації та на відкладеній тестовій вибірці не надто різняця у порівнянні з іншими компонентами.

Для 3 інших моделей результати - невтішні. ## Візуалізація результатів

Зобразмо спочатку загальну картину: гістограму розсіювання з усіма прознозами для всіх моделей:

plot_ly(data = test, x = ~math, y = ~ukr, mode = "markers", type = "scatter", name="EIT results") %>%
  add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 1], name = "pro_ukr",  mode = 'lines+markers') %>%
  add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 2], name = "radical",  mode = 'lines+markers') %>%
  add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 3], name = "opposition",  mode = 'lines+markers') %>%
  add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 4], name = "small",  line=list(color="yellow"), marker=list(color="yellow", line=list(color="yellow")), mode = 'lines+markers') %>%
  add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 5], name = "not_voted",  mode = 'lines+markers') %>%
  layout(title = "Passed EIT results (test dataset)")

Окремо зобразимо для абітурієнтів, які проживають в областях, що ймовірніше наслідують проукраїнську групу:

train <- read.csv("data/computation_results/train2017.csv", row.names = 1, stringsAsFactors = FALSE)
test <- read.csv("data/computation_results/test2017.csv", row.names = 1, stringsAsFactors = FALSE)

prediction <- as.matrix(read.csv("data/computation_results/prediction2017v.csv", row.names = 1, stringsAsFactors = FALSE))
A_coeff <- as.matrix(read.csv("data/computation_results/A_coeff2017v.csv", row.names = 1, stringsAsFactors = FALSE))

Wtest <- as.data.frame(cbind(prediction, test[, -c(1, 2)]))
Wtest[, "max_value"] <- apply(Wtest[, -(1:5)], 1, max)

opposition <- test[Wtest[, "opposition"] == Wtest["max_value"], ]
p1 <- plot_ly(data = opposition, x = ~math, y = ~ukr, mode = "markers", type = "scatter", name = "EIT results in pro_ukr regions") %>%
  add_trace(x = test[rownames(opposition), "math"],
            y = Wtest[rownames(opposition), 3], name = "opposition",  mode = 'markers')
p1

Відмітимо, що основні тренди моделі для проукраїнських регіонів не надто збереглися.

train <- read.csv("data/computation_results/train2017.csv", row.names = 1, stringsAsFactors = FALSE)
test <- read.csv("data/computation_results/test2017.csv", row.names = 1, stringsAsFactors = FALSE)

prediction <- as.matrix(read.csv("data/computation_results/prediction2017v.csv", row.names = 1, stringsAsFactors = FALSE))
A_coeff <- as.matrix(read.csv("data/computation_results/A_coeff2017v.csv", row.names = 1, stringsAsFactors = FALSE))

Wtest <- as.data.frame(cbind(prediction, test[, -c(1, 2)]))
Wtest[, "max_value"] <- apply(Wtest[, -(1:5)], 1, max)

pro_ukr <- test[Wtest[, "pro_ukr"] == Wtest["max_value"], ]
p1 <- plot_ly(data = pro_ukr, x = ~math, y = ~ukr, mode = "markers", type = "scatter", name = "EIT results in pro_ukr regions") %>%
  add_trace(x = test[rownames(pro_ukr), "math"],
            y = Wtest[rownames(pro_ukr), 1], name = "pro_ukr",  mode = 'markers')
p1

train <- read.csv("data/computation_results/train2017.csv", row.names = 1, stringsAsFactors = FALSE)
test <- read.csv("data/computation_results/test2017.csv", row.names = 1, stringsAsFactors = FALSE)

prediction <- as.matrix(read.csv("data/computation_results/prediction2017v.csv", row.names = 1, stringsAsFactors = FALSE))
A_coeff <- as.matrix(read.csv("data/computation_results/A_coeff2017v.csv", row.names = 1, stringsAsFactors = FALSE))

Wtest <- as.data.frame(cbind(prediction, test[, -c(1, 2)]))
Wtest[, "max_value"] <- apply(Wtest[, -(1:5)], 1, max)

not_voted <- test[Wtest[, "not_voted"] == Wtest["max_value"], ]
p1 <- plot_ly(data = not_voted, x = ~math, y = ~ukr, mode = "markers", type = "scatter", name = "EIT results in not_voted regions") %>%
  add_trace(x = test[rownames(not_voted), "math"],
            y = Wtest[rownames(not_voted), 5], name = "not_voted",  mode = 'markers')
p1